home *** CD-ROM | disk | FTP | other *** search
/ SPACE 2 / SPACE - Library 2 - Volume 1.iso / misc / 642 / sdocsata / examples.e < prev    next >
Encoding:
Text File  |  1991-08-06  |  20.2 KB  |  1,126 lines
  1. Schoonschip manual examples.
  2. *end
  3.  
  4. C Example 1
  5.  
  6. Z XX = (a + b)^2
  7.  
  8. *begin
  9.  
  10. B a
  11.  
  12. Z XX= (a + b)^2
  13.  
  14. *end
  15.  
  16. C Example 2.
  17.  
  18. P output   !Only the first two characters of 'output' are significant.
  19.  
  20. Digits 3
  21. Rationalization 3        ! 0 leads to no rationalization attempt.
  22.  
  23. Z XX=3.14 + 1/3*a1 + 7.1357689E20*a2    !E20 or E+20 implies 10^20.
  24.  
  25. *yep
  26.  
  27. Digits            !This restores to the default option, 5 digits.
  28. Rationalization        !Default: 22 digits, with check.
  29.  
  30. *end
  31.  
  32. C Example 3.
  33.  
  34. BLOCK text{xx,yy}
  35. C This argument is substituted in a call: 'xx'.
  36.   This not: xx.
  37.   Arguments can be glued together: 'xx''yy'.
  38. ENDBLOCK
  39.  
  40. BLOCK lines{A,B}
  41. 'A'C line 1
  42. 'A'  line 2
  43. 'B'  line 3
  44. ENDBLOCK
  45.  
  46. C Calling block text:
  47.  
  48. text{yes,sir}
  49.  
  50. DELETE text
  51.  
  52. C Calling lines twice:
  53.  
  54. lines{,~}
  55. lines{~,} 
  56.  
  57. DO var=1,3,2
  58. P input            !Normally only the first DO round is printed.
  59. C value of var: 'var'
  60. ENDDO
  61.  
  62. DO var=-1,0        !The third number is 1 by default.
  63. P input
  64. C value of var: 'var'.
  65. C This may be part of a name, for example XX'var', or an index,
  66.   for example YY(3+'var'). In the last case Schoonschip will interpret
  67.   the % sign as minus if not on a C line and correctly work out the
  68.   number.
  69. ENDDO
  70.  
  71. *end
  72.  
  73. C Example 4.
  74.  
  75. A a1,a2,a3        !Algebraic symbol list.
  76. F f1,f2,f3        !Functions.
  77. I mu,nu            !Indices.
  78. V p,q,k            !Vectors.
  79.  
  80. Z XX = p(mu)*{ a3^-20*a1*q(mu) + a2*a1^7*D(mu,nu)*k(nu)
  81.  
  82.     + a3*f2(mu,k,q)*f1}        !Note the use of {}.
  83.  
  84.     + q(mu)*{a3*a2*a1 + f3*f2*f1}
  85.  
  86. *end
  87.  
  88. C Example 5.
  89.  
  90. A a1,a2=i,a3=c,a4=c,a5=i=3
  91. V p,q=z,k
  92.  
  93. Z XX = Conjg(a1 + a2 + a3 + a4 + a5)
  94. Z YY = Conjg(a3C + a4C)
  95. Z ZZ = (a1 + a5)^5
  96. Z AA = p(mu)*{p(mu) + q(mu)} + D(mu,3)*{p(mu) + q(mu)} + pDk 
  97.  
  98. Oldnew a4C=b4,k=K
  99.  
  100. *end
  101.  
  102. C Example 6.
  103.  
  104. I mu,nu,la,ka
  105.  
  106. Z  XX = D(mu,nu)*f1(a1,a2,mu)*f2(a3,a4,nu)
  107.       + f1(a1,a2,mu)*f2(a3,a4,nu)*D(mu,nu)
  108.  
  109. Z  YY = D(la,ka)*f1(a1,a2,la)*f2(a3,a4,ka)
  110.       + f1(a1,a2,la)*f2(a3,a4,ka)*D(la,ka)
  111.  
  112. Sum,la,ka
  113.  
  114. *end
  115.  
  116. C Example 7.
  117.  
  118. A a1,a2=i,a3=c
  119. F f1,f2=i,f3=c,f4=u
  120.  
  121. Z xx=Conjg{ f1(a1,a2,a3)*f2(a1,a2,a3,f3)*f4(a1,a2,a3) }
  122.  
  123. *end
  124.  
  125. C Example 8.
  126.  
  127. Digits 20
  128. Ratio 0
  129.  
  130. Z xx= DS(J,0,20,(1.),(1/J))
  131.  
  132. *begin
  133.  
  134. Digits
  135. Ratio
  136.  
  137. Z yy = DS(J,0,5,(X^J),(1/J))
  138.  
  139. *end
  140.  
  141. C Example 9.
  142.  
  143. X EXP(y) = DS(J,0,6,(y^J),(1/J))
  144.  
  145. Z COS = { EXP{ (i*x) } + EXP{ -(i*x) } }/2
  146.  
  147. *end
  148.  
  149. C Example 10.
  150.  
  151. B a11,a21,a31
  152.  
  153. D rc(n)=a11,a12,a13,
  154.     a21,a22,a23,
  155.     a31,a32,a33
  156.  
  157. X matrix(n,m)=rc(3*n+m-3)
  158.  
  159. Z xxx = DS{J1,1,3,(DS{J2,1,3,(DS{J3,1,3,
  160.     { DP(J1,J2,J3)*matrix(J1,1)*matrix(J2,2)*matrix(J3,3) }
  161.     } ) } ) }
  162.  
  163. *begin
  164.  
  165. B a11,a21,a31,b11,b21,b31
  166.  
  167. D ra(n)=a11,a12,a13,
  168.     a21,a22,a23,
  169.     a31,a32,a33
  170.  
  171. D rb(n)=b11,b12,b13,
  172.     b21,b22,b23,
  173.     b31,b32,b33
  174.  
  175. X matrix(n,m,ra)=ra(3*n+m-3)
  176.  
  177. X DET(ra) = DS{J1,1,3,(DS{J2,1,3,(DS{J3,1,3,
  178.     { DP(J1,J2,J3)*matrix(J1,1,ra)*matrix(J2,2,ra)*matrix(J3,3,ra) }
  179.     } ) } ) }
  180.  
  181. Z detb=DET(rb)
  182.  
  183. Z deta=DET(ra)
  184.  
  185. *begin
  186.  
  187. F f1
  188.  
  189. B b11,b21,b31
  190.  
  191. D ra(n)=(a1-a2),0,
  192.     0,(a1+a2)
  193.  
  194. D rb(n)=b11,b12,b13,
  195.     b21,b22,b23,
  196.     b31,b32,b33
  197.  
  198. D rc(n)=c11
  199. X matrix(k,n,m,f1)=f1(k*n+m-k)
  200.  
  201. D DET(n,f1) = f1(1),
  202.    DS{J1,1,2,(DS{J2,1,2,
  203.      { DP(J1,J2)*matrix(n,J1,1,f1)*matrix(n,J2,2,f1) }
  204.      } ) }
  205.   ,
  206.    DS{K1,1,3,(DS{K2,1,3,(DS{K3,1,3,
  207.      { DP(K1,K2,K3)*matrix(n,K1,1,f1)*matrix(n,K2,2,f1)*matrix(n,K3,3,f1) }
  208.     } ) } ) }
  209.  
  210. Z deta=DET(2,ra)
  211.  
  212. Z detb=DET(3,rb)
  213.  
  214. Z detc=DET(1,rc)
  215.  
  216. *end
  217.  
  218. C Example 11.
  219.  
  220. T List(n)=a,b,c
  221.  
  222. F f
  223.  
  224. Z compl = DS(J,1,3,{ f(List(J)) } )
  225.  
  226. *begin
  227.  
  228. A x=c
  229.  
  230. T r1(n)=a11,a12,a13
  231. T r2(n)=a21,a22,a23
  232. T r3(n)=a31,a32,a33
  233.  
  234. T matr(n,m)=r1,r2,r3
  235.  
  236. T weird(n,a1,a2)=Conjg(a1+a2),Integ(3*a2)
  237.  
  238. X XX(a1,a2)=a2*a1
  239.  
  240. Z sqa13 = DS(J,1,3,{f1(matr(1,J))*f1(matr(J,3))} )
  241.  
  242. Z weirdo=XX(weird(1,x,7),weird(2,x,7))
  243.  
  244. *end
  245.  
  246. C Example 12.
  247.  
  248. P lists
  249.  
  250. Common yyy,ccc(0)
  251.  
  252. Z xxx(a,b)=(a+b)^2
  253. Z ccc(3,a,b)=(a+b)^3
  254. Z ccc(4,a,b)=(a+b)^4
  255.  
  256. Keep xxx
  257.  
  258. *next
  259.  
  260. F f1,f2
  261. V p
  262.  
  263. Z yyy(p)=xxx(c,d) + p(nu)*f2(mu)
  264.  
  265. Sum mu,nu
  266.  
  267. *begin
  268. Write fileC
  269. *begin
  270.  
  271. V q
  272.  
  273. Z zzz=a1*yyy(q)
  274.  
  275. Z Abc=xxx(a,b)
  276.  
  277. Delete ccc,yyy        !The actual delete when this section is done.
  278.  
  279. *begin
  280.  
  281. Z xyz=ccc(3,e,f)
  282.  
  283. *end        !The reader may want to delete fileC at this point.
  284.  
  285. C Example 13.
  286.  
  287. Z integr = a*x^2 +b*x
  288.  
  289. Id,x^n~ = x2^(n+1)/(n+1) - x1^(n+1)/(n+1)
  290.  
  291. *begin
  292.  
  293. Z integr = (a*x^2 + b*x + c)*dx
  294.  
  295. Id,x^n~*dx = x2^(n+1)/(n+1) - x1^(n+1)/(n+1)
  296.  
  297. Al,dx = x2 - x1
  298.  
  299. *begin
  300.  
  301. C Method to integrate expressions of the form x^n*sin(x) or x^n*cos(x).
  302. C The method is based on the equation (n even):
  303. C Integral( x^n*sin(x)) = - n! * cos(x)
  304.      + Integral(  sin(x) * n * { x^(n-1) - (n-2)*(n-1)*x^(n-3) + ...} )
  305.      + Integral( -cos(x) * { x^n - (n-1)*n*x^(n-2) + ...} )
  306.   and similar equations for odd n and also cos(x).
  307.   The X-expressions Coefa and Coefb corresponds to the above sums enclosed
  308.   in curly brackets.
  309.   Below the function DK is used to separate even and odd n cases. Remember
  310.   that Integ converts and truncates its argument to a short integer in the
  311.   range -128, 127.
  312.  
  313. X Coefa(n,m,x)=DS(j,0,m,{ x^(n-2*j-1) } , { -(n-2*j)*(n-2*j+1) } )
  314. X Coefb(n,m,x)=DS(k,0,m,{ x^(n-2*k) } , { -(n-2*k+2)*(n-2*k+1) } )
  315.  
  316. X Icos(n,x) =
  317.    DK(2*Integ(n/2),n)*{
  318.     - DB(n)*sin(x)                ! n even.
  319.     + cos(x)*n*Coefa(n,(n-1)/2,x)
  320.     + sin(x)*Coefb(n,(n-1)/2-1,x)  }
  321.  
  322.  + DK(2*Integ(n/2),n+1)*{
  323.       DB(n)*cos(x)                ! n odd.
  324.     + cos(x)*n*Coefa(n,(n-1)/2-1,x)
  325.     + sin(x)*Coefb(n,(n-1)/2,x)    }
  326.  
  327. X Isin(n,x) =
  328.    DK(2*Integ(n/2),n)*{
  329.       DB(n)*cos(x)                ! n even.
  330.     + sin(x)*n*Coefa(n,(n-1)/2,x)
  331.     - cos(x)*Coefb(n,(n-1)/2-1,x)  }
  332.  
  333.  + DK(2*Integ(n/2),n+1)*{
  334.     + DB(n)*sin(x)                ! n odd.
  335.     + sin(x)*n*Coefa(n,(n-1)/2-1,x)
  336.     - cos(x)*Coefb(n,(n-1)/2,x)    }
  337.  
  338. Z xx = x^8*cos(x)
  339.  
  340. Id,x^n~*cos(x) = Icos(n,x2) - Icos(n,x1)
  341. Al,x^n~*sin(x) = Isin(n,x2) - Isin(n,x1)
  342.  
  343. *end
  344.  
  345. C Example 14.
  346.  
  347. V p,q
  348. A a,b,c,d,e,f,g,h
  349. I mu,nu
  350. F f1,f2,f3
  351.  
  352. X expr=f1(a,b,p)*f2(a,c,q)*f3(d,e)*f1(g,h)*
  353.     { a^7 + a^-7 + a^2 + pDq^2 + pDp + p(mu) + p(nu) }
  354.  
  355. *fix
  356.  
  357. C Class 0, no keyword.
  358.  
  359. Z xxx=expr
  360.  
  361. Id,f1(a1~,a2~)=a1^10*a2^20
  362.  
  363. *begin
  364.  
  365. C Class 0, keyword Always.
  366.  
  367. Z xxx=expr
  368.  
  369. Id,Always,f1(a1~,a2~)=a1^10*a2^20
  370.  
  371. *begin
  372.  
  373. C Class 1, keyword Multi.
  374.  
  375. Z xxx=expr
  376.  
  377. Id,Multi,a^3 = xyz + hij
  378.  
  379. *begin
  380.  
  381. C Class 2, exponent 1, no keyword.
  382.  
  383. Z xxx=expr
  384.  
  385. Id,pDq = XYZ + HIJ
  386.  
  387. *begin
  388.  
  389. C Class 3, no keyword.
  390.  
  391. Z xxx=expr
  392.  
  393. Id,a^2 = a1^7/15
  394.  
  395. *begin
  396.  
  397. C Class 5.
  398.  
  399. Z xxx=expr
  400.  
  401. Id,p(mu) = - q(mu)
  402.  
  403. *begin
  404.  
  405. C Class 6.
  406.  
  407. Z xxx=expr
  408.  
  409. Id,p(mu~) = - q(mu)
  410.  
  411. *begin
  412.  
  413. C Class 7, keyword Funct.
  414.  
  415. Z xxx=expr
  416.  
  417. Id,Funct,a = a27
  418.  
  419. Id,f1(a1~,a2~,a3~) = 200*a1*a2*a3
  420. Al,f2(a1~,a2~,a3~) = a1^10*a2^11*a3^12
  421.  
  422. *begin
  423.  
  424. C Class 8, keyword Once.
  425.  
  426. Z xxx=expr
  427.  
  428. Id,Once,a^2 = XXX
  429.  
  430. *begin
  431.  
  432. C Class 9.
  433.  
  434. Z xxx=expr
  435.  
  436. Id,a^2*f2(a1~,c,p~) = F2(a1,c,p)
  437.  
  438. *begin
  439.  
  440. C Class 10.
  441.  
  442. Z xxx=expr
  443.  
  444. Id,f1~(a,b~,p~) = F(a,b,p)
  445.  
  446. *begin
  447.  
  448. C Class 11, keyword Adiso.
  449.  
  450. Z xxx = expr
  451.  
  452. Id,Adiso,f1(g,h)*f1(a~,b~,c~) = F1(a,b,c,g,h)
  453.  
  454. *begin
  455.  
  456. C Class 12, keyword Ainbe.
  457.  
  458. Z xxx=expr
  459.  
  460. Id,Ainbe,f1(g,h)*f1(a~,b~,c~) = F1(a,b,c,g,h)
  461. Id,Ainbe,f1(a~,b~,c~)*f1(g,h) = F2(a,b,c,g,h) 
  462.  
  463. *begin
  464.  
  465. C Class 13.
  466.  
  467. Z xxx=expr
  468.  
  469. Id,f1(x~,b,p)*f2(x~,c,q) = F(x,b,p,c,q)
  470.  
  471. *begin
  472.  
  473. C Class 15, keyword Dotpr.
  474.  
  475. Z xxx=expr
  476.  
  477. Id,Dotpr,p(mu~) = - q(mu)
  478.  
  479. *begin
  480.  
  481. C Class 16, keyword Funct.
  482.  
  483. Z xxx=expr
  484.  
  485. Id,Funct,p(mu~) = - F1(a,b,mu)
  486.  
  487. *end
  488.  
  489. C Example 15.
  490.  
  491.   Integration of polynomium.
  492.  
  493. A a,b,c,d,e,x
  494.  
  495. Z xxx = a*x^2 + b*x + c + d/x + e/x^2
  496.  
  497. IF NOT x^-1=[Log(x)]
  498. AND NOT x^n~=x^(n+1)/(n+1)
  499. Al,Addfa,x
  500. ENDIF
  501.  
  502. *end
  503.  
  504. C Example 16.
  505.  
  506. V p,q
  507. F F
  508. I mu=N,nu=N
  509.  
  510. Z xx= pDp^2 * pDq^3 * F1(p,p,q,p) * p(mu) * q(nu)
  511.  
  512. Id,All,p,N,F
  513. P output
  514. *yep
  515. C Showing the dimensionality of the created indices:
  516.  
  517. Id,F(i1~,i2~,i3~,i4~,i5~,i6~,i7~,i8~,i9~,i10~,i11~)=
  518.  F(i1,i2,i3,i6,i7,i8,i9,i10,i11)*D(i4,i5)
  519. *begin
  520.  
  521. C A more realistic example.
  522.  
  523. I al,be,mu,nu
  524. A N,N_
  525. V p,k
  526. F F,F20,F22
  527.  
  528. Z xx = G(1,1,al,be,p,al,be,p)
  529.  
  530. Id,All,p,N,F
  531. P output
  532. *yep
  533. Id,F(i1~,i2~) = D(i1,i2)*F20(k) + k(i1)*k(i2)*F22(k)
  534. *end
  535.  
  536. C Example 17.
  537.  
  538. V p,q,k
  539. I mu
  540. A a,b,c,d,e,f,g,h
  541. F F1,F2,F3
  542.  
  543. Z xx = F1(e,d,c,b,a) + F2(e,d,c,b,a) + F3(e,d,c,b,a)
  544.     + F2(-125,-30,-1,0,30,125)
  545.     + F2(pDq,pDk,kDq,p(3),q(2),F2,7,mu,p,a)
  546.     + f1(e,f,g,d,a,b,c,h,k)
  547.  
  548. Id,Asymm,F1,2,3,4,F2,F3,4,5
  549. Al,Asymm,f1:1-3:5-7:
  550.  
  551. *end
  552.  
  553. C Example 18.
  554.  
  555. A a1,a2,a3,a4,alt1,agt2,aeq15,ai12
  556.  
  557. Z xxx= 0.5*a1 + 1.5*a2 + 2.5*a3 + 1.E-20*a4
  558.  
  559. IF Coef,lt,1.
  560. Al,Addfa,alt1
  561. ENDIF
  562.  
  563. IF Coef,gt,2.
  564. Al,Addfa,agt2
  565. ENDIF
  566.  
  567. IF Coef,eq,1.5
  568. Al,Addfa,aeq15
  569. ENDIF
  570.  
  571. IF Coef,ib,1,2
  572. Al,Addfa,ai12
  573. ENDIF
  574.  
  575. IF Coef,lt,1E-15    !Delete the term if coefficient < 1E-15.
  576. Al,Addfa,0
  577. ENDIF
  578.  
  579. *end
  580.  
  581. C Example 19.
  582.  
  583. A a1,a2,a3,a4,a5
  584. F F1,F2,F3
  585.  
  586. Z xx= F3(a1,a2)*F1(a1,a2,a3)*F1(a4,a5)
  587.  + F1(a1,a2,a4)*F2(a1,a2,a3)*F1(a1,a2,a3)
  588.  + F3(a1,a4)*F2(a1,a2,a3)*F1(a1,a2,a4)*F2(a3,a4,a5)
  589.  
  590. Id,Commu,F1,F3
  591.  
  592. *end
  593.  
  594. C Example 20.
  595.  
  596. F F1,F2
  597. A a,b,c,d
  598.  
  599. Z xx =      F1("c,"b,"a,/,"e,*,a1)
  600.     + F1("c,"b,"a,/,"e,*,a1,*,xy)
  601.     + F1("c,"b,"a,/,"e,*,a1,*,xy,*,a2)
  602.     + F1("c,"b,"a,/,"e,*,a1,*,xy,*,a2,*,xz)
  603.  
  604.     + F2("a,"b,/,"F,*,x,*,y,*,z1)
  605.     + F2("b,"a,/,"F,*,x,*,y,*,z1)
  606.     + F2("F,"c,"b,"a,*,*,z,*,y,*,x)
  607.     + F2("F,"c,"b,"a,*,*,z,*,yp,*,x)
  608.  
  609.  
  610. Id,Compo,<AVIFXA>,F1,<FF>,F2
  611. Id,Always,F2(F1~,a~,b~,c~,d~,e~) = F1(a,b,c)
  612. Id,F1~(a~,b~,z1) = F1(a,b)
  613. *begin
  614. C
  615.            j1
  616.     i1   ---0-------0---  i3
  617.          j4    |    |j2
  618.         |    |
  619.     i2   ---0-------0---  i4
  620.            j3
  621.  
  622.  
  623. T Ch(n)="A,"B,"C
  624. T Cg(n)="A,"C,"B
  625.  
  626. F BC,ABC
  627. A i1,i2,i3,i4,j1,j2,j3,j4
  628.  
  629. Z solu = square("A,"A,"A,"A)
  630.  
  631. Id,square(c1~,c2~,c3~,c4~) = DS{L1,1,3,(DS{L2,1,3,
  632.     (DS{L3,1,3,(DS{L4,1,3,(
  633.         v3(c1,Ch(L1),Cg(L4),*,i1,*,-j1,*,j4)*
  634.         con(Ch(L1),Cg(L1),*,j1,*,-j1)*
  635.         v3(c3,Cg(L1),Ch(L2),*,-i3,*,j1,*,-j2)*
  636.         con(Ch(L2),Cg(L2),*,j2,*,-j2)*
  637.         v3(c4,Cg(L2),Ch(L3),*,-i4,*,j2,*,-j3)*
  638.         con(Ch(L3),Cg(L3),*,j3,*,-j3)*
  639.         v3(c2,Cg(L3),Ch(L4),*,i2,*,j3,*,-j4)*
  640.         con(Ch(L4),Cg(L4),*,j4,*,-j4)
  641.         ) } ) } ) } ) }
  642.  
  643. Id,Compo,<F>,v3,con
  644. Id,v3(f1~,a1~,a2~,a3~) = f1(a1,a2,a3)
  645. Al,con(f1~,a1~,a2~) = f1(a1)
  646.  
  647. *end
  648.  
  649. C Example 20.
  650.  
  651. A x,y
  652. T tt(n) = "1,"2,"3,"4,"5
  653.  
  654. C Here a complicated way to make
  655. C    pow(x) = a1*x + a2*x^2 + a3*x^3 + a4*x^4 + a5*x^5.
  656.  
  657. Z pow(x) = DS(J,1,5,{f1(/,"a,tt(J))*x^J})
  658.  
  659. Id,Compo,f1
  660. Id,f1(y~) = y
  661.  
  662. Keep pow
  663. *next
  664.  
  665. Z xx = pow(x)
  666.  
  667. Id,Count,3,x,1
  668.  
  669. Keep pow
  670. *next
  671.  
  672. Z xx = pow(y)
  673.  
  674. Id,Count,f1,y,2,a3,10
  675.  
  676. Keep pow
  677. *next
  678. V p,q
  679. A AA
  680.  
  681. Z xx = pow(y)*f1(a2,a3)*pDq^2
  682.  
  683. Id,Count,AA,y,2,f1,-4,p,1,q,3
  684.  
  685. *end
  686.  
  687. C Example 22.
  688.  
  689. A a,b,c,d,e
  690.  
  691. Z xxx = F1(e,d,c,b,a) + F2(e,d,c,b,a) +F3(e,d,c,b,a)
  692.  
  693. Id,Cyclic,F1,2,5,4
  694. Id,Symme,F2,F3,2,3,4
  695.  
  696. *end
  697.  
  698. C Example 23.
  699.  
  700. V p,q,k,pp,qp,kp
  701.  
  702. Z xxx = Epf(k,p,q)*Epf(kp,pp,qp) + Epf(i1,i2,i3,i4)*Epf(i1,i2,j3,j4)
  703.  
  704. Id,Epfred
  705.  
  706. *end
  707.  
  708. C Example 24.
  709.  
  710. F f1,f2,f3,f4
  711.  
  712. Z xxx = f1(-a1,-a2,a3,-a4) + f2(-a1,-a2,a3,-a4) + f3(-a1,-a2,a3,-a4)
  713.     + f4(-a1,-a2,a3,-a4)
  714.  
  715. Id,Even,f1,2,3,f2
  716. Id,Odd,f3,2,3,f4
  717.  
  718. *end
  719.  
  720. C Example 25.
  721.  
  722. V p,q
  723.  
  724. Z xxx = a1^2*pDq^3/p(4)
  725.  
  726. Id,Numer,a1,2,pDq,1.E10,p(4),1.E5
  727.  
  728. *end
  729.  
  730. C Example 26.
  731.  
  732. A a1,a2,a2ma1
  733. F f1
  734. B b2,b3,b4,b5
  735.  
  736. Z xx=f1(8,4)
  737.  
  738. Id,f1(n~,m~)=
  739.     { b2*a1^n*a2^m 
  740.     + b3*a1^-n*a2^m 
  741.     + b4*a1^n*a2^-m
  742.     + b5*a1^-n*a2^-m
  743.     }
  744.  
  745. Id,Ratio,a1,a2,a2ma1
  746. *begin
  747.  
  748. C Here the use of Ratio with [] names.
  749.  
  750. B [b-a]
  751.  
  752. Z xxx = 1/[x+a]^3 * 1/[x+b]^2
  753.  
  754. Id,Ratio,[x+a],[x+b],[b-a]
  755. P output
  756. *yep
  757. C Just checking...
  758.  
  759. Id,[x+a]^n~ = (x+a)^(3+n)*(x+b)^2/[x+a]^3/[x+b]^2
  760. Al,[x+b]^n~ = (x+a)^3*(x+b)^(2+n)/[x+a]^3/[x+b]^2
  761. Id,b = [b-a] + a
  762. *end
  763.  
  764. C Example 27. Muon decay.
  765.  
  766.   Masses and momenta:    muon Mm,k; electron Me,p;
  767.             anti-e-neutrino qp; mu-neutrino q. 
  768.   The length of the 3-dimensional parts of p, q and qp are
  769.   denoted by pl, ql and qpl. The variable z is the cosine
  770.   of the angle between p and q (3-dimensional parts).
  771.   Also k0, p0 etc are the energies of the particles. Note
  772.   that k4=i*k0 etc. Note also that q0=ql, qp0=qpl since
  773.   the neutrino masses are taken to be zero.
  774.  
  775. C Further quantities: Pi=3.14... and Mpr is the proton mass.
  776.  
  777. C The evaluation is in the muon rest-system, thus kl=0 and
  778.   k0=Mm.
  779.  
  780. A Me,Mm,ql,qpl,pl,p0,Pi,b,z,bp,Mpr,Al
  781. V p,q,qp,k
  782. I mu,mup,i1,i2
  783. F Dia=u
  784.  
  785. Z     Rate = Dia(mu,1)*Conjg(Dia(mup,10))/2
  786.  
  787. C Factor 1/2 for averaging over muon polarizations.
  788.  
  789. Id,Dia(mu~,n~) = Ubg(n,Me,p)*G(n,n+1,mu,G6)*Ug(n+1,0,qp)
  790.              *Ubg(n+2,0,q)*G(n+2,n+3,mu,G6)*Ug(n+3,Mm,k)
  791. Id,Spin,p,q,k,qp
  792.  
  793. C There are two separate traces here. They may be unified to
  794.   advantage as they have two indices in common (mu and mup).
  795.  
  796. Id,Gammas,"U
  797.  
  798. C Apply conservation of four-momentum  qp=k-p-q  and
  799.   the mass-shell relations kDk=-Mm^2, pDp=-Me^2 and
  800.   qDq=qpDqp=0.
  801.  
  802. Id,qpDq=kDq-pDq
  803. Al,qpDk=-Mm^2-pDk-qDk
  804. Al,qpDp=pDk+Me^2-qDp
  805. Id,Funct,qp(mu~) = k(mu) - p(mu) - q(mu)
  806.  
  807. *yep
  808.  
  809. Id,pDq=pl*ql*z-p0*ql
  810. Al,pDk=-Mm*p0
  811. Al,qDk=-ql*Mm
  812. Al,kDk=-Mm^2
  813. Al,pDp=-Me^2
  814.  
  815. C Integration over all angles of the vector p
  816.   gives a factor 4*Pi:
  817.  
  818. Id,Addfa,4*Pi
  819.  
  820. C The integration over the azimuthal angle of q
  821.   gives a factor 2*Pi. Including the various factors
  822.   1/(2k0)..., a factor (2*Pi)^4/{(2*Pi)^3)}^3
  823.   and the factors from going to polar coordinates
  824.   for q and p gives:
  825.  
  826. Id,Addfa,2*Pi/16/Mm/p0/ql/qpl*pl^2*ql^2/32/Pi^5
  827.  
  828. C The remaining delta function is delta(k0-p0-q0-qp0).
  829.   For given p0=Sqrt(pl^2+Me^2) and q0=ql this may be
  830.   solved for z using qp0=qpl=length(p+q)=
  831.   Sqrt(pl^2+ql^2+2*pl*ql*z) where z is the cosine of
  832.   the angle of q with the third axis, taken along p. 
  833.   In this last comment p and q are three-vectors.
  834.  
  835. C Note that integration over the delta-function gives
  836.   a factor 1/Abs(F), where F is the derivative of the
  837.   argument of the delta function with respect to z.
  838.   This F is equal to ql*pl/qpl.
  839.  
  840. Id,z=(0.5*Me^2+0.5*Mm^2-Mm*p0-ql*Mm+ql*p0)/pl/ql
  841. Al,Addfa,qpl/ql/pl
  842.  
  843. *yep
  844.  
  845. B Pi
  846.  
  847. C The next integration is over the length ql. The
  848.   endvalues are denoted by qma and qmi, with
  849.   qma=(Mm-p0+pl)/2 and qmi=(Mm-p0-pl)/2. These
  850.   values obtain for configurations where the electron
  851.   and both neutrinos are aligned. The maximum value
  852.   obtains if the momentum qp is in the direction
  853.   of the electron momentum, while q is pointed in the
  854.   opposite direction. Then pl+qpl=ql and q0=Mm-p0-qp0,
  855.   and one solves ql(=q0)=(Mm-p0+pl)/2.
  856.  
  857. C Actually, ql^-1 is not occurring here, but for
  858.   completeness we include it in the integration.
  859.  
  860. IF NOT ql^-1=[Log(qa/qb)]
  861. AND NOT ql^n~=qma^(n+1)/(n+1)-qmi^(n+1)/(n+1)
  862. Al,Addfa,qma-qmi
  863. ENDIF
  864.  
  865. Id,qma=0.5*Mm-0.5*p0+0.5*pl
  866. Al,qmi=0.5*Mm-0.5*p0-0.5*pl
  867.  
  868. *yep
  869.  
  870. Id,Multi,p0^2 = pl^2+Me^2
  871. Id,p0 = (pl^2+Me^2)/p0
  872.  
  873. C The answer contains essentially one kinematic variable,
  874.   namely pl (or p0). The spectrum that follows from this
  875.   equation is typical for a V-A theory. Using a variable
  876.   x (shown below) this spectrum follows an equation derived
  877.   by Michel provided a parameter ro in that equation, called
  878.   Michel parameter, is set to 3/4.
  879.  
  880. P output
  881. *yep
  882.  
  883. C The final integration is over pl. A new variable
  884.   x defined by x=pl+p0, with p0=Sqrt(pl^2+Me^2)
  885.   is used in the following. One has dx/d(pl)=1+p/p0
  886.   and this may be rewritten as d(pl)/p0 = dx/x.
  887.   Note that pl=(x+Me^2/x)/2.
  888.   There are terms with and without p0^-1, and they
  889.   are treated separately.
  890.   The variable plmx is the maximum value of pl, i.e.
  891.   plmx=(Mm-Me^2/Mm)/2. This is achieved if both neutrinos
  892.   are aligned and opposite to the electron. The
  893.   energy of the neutrinos is then simply pl, since
  894.   the sum of the momenta must be equal to the electron
  895.   momentum. In that case therefore Mm=p0+pl, which gives
  896.   the stated value for plmx. Note that Mm is the endpoint
  897.   value for x, corresponding to maximum pl.
  898.   The minimum value for pl is 0, for x it is Me.
  899.  
  900. IF p0^-1=1/x
  901. Al,pl^n~=(x/2-Me^2/x/2)^n
  902. Id,x^-1=[Log(Mm/Me)]
  903. Al,x^n~=Mm^(n+1)/(n+1)-Me^(n+1)/(n+1)
  904. ELSE
  905. ..IF NOT pl^n~=plmx^(n+1)/(n+1)
  906. ..Al,Addfa,plmx
  907. ..ENDIF
  908. ENDIF
  909.  
  910. Id,plmx=Mm/2-Me^2/Mm/2
  911.  
  912. P output
  913. *yep
  914.  
  915. Digits 7
  916.  
  917. C Numerical evaluation. One usually introduces the coupling
  918.   constant g/Sqrt(2)/Mpr^2, and then g is to be taken such
  919.   that the muon lifetime comes out to be 2.197134E-6.
  920.   Note also the factor (h/(2*Pi))^-1 to convert Mev to 1/seconds.
  921.  
  922. Id,Addfa,g^2/6.582173E-22/2/Mpr^4
  923. Id,Numer,Mm,105.65946,Me,0.5110034,[Log(Mm/Me)],5.33160,Mpr,938.2796
  924. Al,Numer,Pi,3.141592653589793238
  925.  
  926. C Multiply the Rate (= inverse lifetime) with the experimental
  927.   lifetime. Then g can be solved from the requirement that
  928.   the result must be one. This gives g = 1.024E-5, which
  929.   corresponds to the well-known result for the fermi-coupling
  930.   constant:
  931.     Gf = 1.02 * 10^-5 /Sqrt(2) /Mpr^2 .
  932.  
  933. Id,Addfa,2.197134E-6
  934.  
  935. P output
  936. *yep
  937.  
  938. Id,g = 1.0246275E-5
  939.  
  940. *end
  941.  
  942. C Example 28.
  943.  
  944. C The problem is to make a series expansion of a
  945.   complicated expression in terms of quantities a and b.
  946.   The expression was given in terms of the quantities
  947.   called X1,X2,Omx1a,X2b,Bmx2a,Fac1 and Fac2 here below.
  948.  
  949. C The problem was solved in two separate runs. Here the
  950.   first run is shown. The object is to give the series
  951.   expansion for the quantities mentioned. Up to order 7
  952.   in a and/or b is required. The quantity ep is an expansion
  953.   parameter, essentially counting the order in a and b.
  954.   Terms ep^8 and higher are to be ignored, the sooner the
  955.   better. For this the assignment ep=8 in the A-list is
  956.   the best tool.
  957.  
  958. C The calculation proceeds by building up one quantity
  959.   after the next, step by step.
  960.  
  961. Common X1,X2,Omx1a,X2b,Bmx2a,Fac1,Fac2
  962.  
  963. A ep=8,a,b
  964. X X(a)=ep^2*(a^2+b^2-2*a*b)-2*ep*(a+b)
  965.  
  966. *fix
  967.  
  968. C First certain roots called WO and WOI are needed.
  969.   These are of the form Sqrt(1+XX) and 1/Sqrt(1+XX).
  970.   The quantity XX itself is a function of a and b,
  971.   given above as X(a). In the following DO-loop the
  972.   powers of X(a) are computed:
  973.   XX(1)=X(a), XX(2)=X(a)^2 etc.
  974.  
  975. Z XX(0)=1.
  976.  
  977. Keep XX
  978. *next
  979.  
  980. DO K1=1,7
  981. B ep
  982. Z XX('K1')=X(a)*XX('K1'-1)
  983. Nprint XX('K1')
  984. Keep XX
  985. *next
  986. ENDDO
  987.  
  988. C Now the roots WO and WOI can be computed. The expansions
  989.   used here are simply the expansions for Sqrt(1+x) and
  990.   1/Sqrt(1+x).
  991.  
  992. B ep
  993.  
  994. Z WO = 1 + XX(1)/2 - XX(2)/8 + XX(3)/16 - 5*XX(4)/128
  995.     + 7*XX(5)/256 - 21*XX(6)/1024 + 33/2048*XX(7)
  996. Z WOI = 1 - XX(1)/2 + 3*XX(2)/8 - 5*XX(3)/16 + 35*XX(4)/128
  997.     - 63*XX(5)/256 + 231*XX(6)/1024 - 429/2048*XX(7)
  998.  
  999. Keep WO,WOI
  1000. *next
  1001.  
  1002. B ep
  1003.  
  1004. Z X1=0.5+0.5*b*ep-0.5*a*ep+0.5*WO
  1005. Z X2=0.5+0.5*b*ep-0.5*a*ep-0.5*WO
  1006.  
  1007. Keep WOI
  1008. *next
  1009.  
  1010. B ep
  1011.  
  1012. Z Omx1a=1/a/ep-X1/a/ep
  1013. Z X2b=X2/b/ep
  1014. Z Bmx2a=b/a-X2/a/ep
  1015.  
  1016. Keep WOI
  1017. *next
  1018.  
  1019. B ep
  1020. A ep=7
  1021.  
  1022. C For Fac1 and Fac2 only up to ep^6 is needed.
  1023.  
  1024. Z Fac1=(b*ep-X1)*WOI
  1025. Z Fac2=(X2-b*ep)*WOI
  1026.  
  1027. *begin
  1028.  
  1029. C Wrinting of common files. The printed list shows which
  1030.   common files were in existence here; a common file
  1031.   has =C after its name. The X-expression X uses two
  1032.   levels and has =X2 after its name. It is not written.
  1033.  
  1034. P lists
  1035. Write expP
  1036. *end
  1037.  
  1038. C Example 29.
  1039.  
  1040. C In this second part (Example 28 is first part)
  1041.   we will not show the whole calculation of all
  1042.   expressions as they were needed, but just
  1043.   the evaluation of the first one, which is
  1044.   enough to clarify the procedure. First the
  1045.   output of example 27 must be entered.
  1046.  
  1047. Enter expP
  1048. A ep=7,a,b
  1049. *fix
  1050.  
  1051. B ep
  1052.  
  1053. Z Arg(1) = X2
  1054. Z Arg(3) = (b-a)*X2/b
  1055. Z Arg(4) = 1 - X2b
  1056. Z Arg(5) = Omx1a*X2b - 1
  1057. Z Arg(6) = (b-a)*(1-X1)/a
  1058. Z Arg(7) = (a-b)*ep*Omx1a
  1059. Z Arg(8) = Bmx2a
  1060.  
  1061. Keep Arg
  1062. *next
  1063.  
  1064. B ep
  1065.  
  1066. C The Sp are certain functions of a variable x.
  1067.   Sp(1) = Log(1-x)/x, higher Sp are obtained by
  1068.   iteration, involving differentiation.
  1069.   In the end x must be set to (a-b)/a.
  1070.   The XX2 are powers of -Arg(1).
  1071.  
  1072. Z Sp(1) = - Lomx/x
  1073. Z XX2(1) = - Arg(1)
  1074.  
  1075. Nprint Sp,XX2
  1076. Keep Sp,XX2,Arg
  1077. *next
  1078.  
  1079. B ep
  1080.  
  1081. DO K1=2,6
  1082. Z Sp('K1') = Difx*Sp('K1'-1)/'K1'
  1083. Z XX2('K1') = - XX2('K1'-1)*Arg(1)
  1084.  
  1085. Id,Difx*Lomx = - 1/Omx + Lomx*Difx
  1086. Id,Difx*Omx^n~ = - n*Omx^(n-1) + Omx^n*Difx
  1087. Id,Difx*x^n~ = n*x^(n-1) + x^n*Difx
  1088. Id,Difx = 0
  1089.  
  1090. C Note that only negative exponents of x occur here.
  1091. Id,x^n~*Omx^-1 = DS(J1,1,-n,(x^-J1)) + 1/Omx
  1092.  
  1093. Nprint Sp,XX2
  1094. Keep Sp,XX2,Arg
  1095. *next
  1096. ENDDO
  1097.  
  1098. B Lga,Lgb
  1099.  
  1100. DO K2=1,6
  1101. Z Sp('K2') = Sp('K2')*x^'K2'
  1102. ENDDO
  1103.  
  1104. C Omx = 1-x appears only with negative exponent. It
  1105.   is rewritten in terms of x/(1-x) = Xomx, which
  1106.   is possible here since the exponent of x is here
  1107.   always larger than minus the exponent of Omx.
  1108.  
  1109. Id,x^n~*Omx^m~ = Xomx^-m*x^(n+m)
  1110. Id,x = (a-b)/a
  1111. Al,Xomx = (a-b)/b
  1112. Id,Lomx = Lgb - Lga
  1113.  
  1114. Nprint Sp
  1115. Keep Sp,XX2,Arg
  1116. *next
  1117.  
  1118. C Now finally the first of the desired expressions
  1119.   is worked out.
  1120.  
  1121. B ep,Lga,Lgb
  1122.  
  1123. Z Exp(1) = DS(J3,1,6,{Sp(J3)*XX2(J3)})
  1124.  
  1125. *end
  1126. ə